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LONG-TERM  GOALS 

The  long  term  goal  is  to  develop  a  robust  multi-phase,  multi-class  numerical  modeling  framework  for 
both  cohesive  and  non-cohesive  sediment  transport  in  the  fluvial  and  coastal  environments. 

OBJECTIVES 

The  objectives  of  the  present  study  focus  on  extending  a  fluid  mud  model  for  boundary-layer-driven 
and  gravity-driven  transport  (Hsu  et  al.  2007)  with  several  new  capabilities.  Specifically,  the  objectives 
are  to: 

•  Carry  out  model  development  and  model-data  comparisons  with  laboratory/field  observed  fluid 
mud  processes,  specifically  in  conjunction  with  new  MURI  initiative  for  understanding  wave- 
mud  interaction. 

•  Extend  the  existing  fluid  mud  model  with  a  flocculation  module  and  a  bed  module  to  enable 
direct  modeling  on  consolidation/fluidization.  Extend  the  model  to  multi-dimensional  for 
various  coastal  and  estuarine  applications. 

•  In  conjunction  to  NOPP-Community  Sediment  Transport  Model  initiative,  the  fluid  model  is 
utilized  to  provide  parameterizations  for  wave-boundary-layer-scale  transport  processes. 

APPROACH 

It  is  important  to  understand  the  fate  of  terrestrial  sediment  into  the  coastal  ocean,  because  it 
determines  for  example  the  seabed  properties  and  the  turbidity  of  the  water  column.  Several  recent 
initiatives,  such  as  NOPP  Community  Sediment  Transport  Model,  MURI-  Understading  Wave-Mud 
Interaction  and  Tidal  Flats  DRI  have  put  forward  new  outstanding  science  and  technical  questions  for 
coastal  sediment.  The  success  of  these  new  studies  inevitably  depends  on  our  level  of  understandings 
on  the  crucial  small-scale  mechanisms  that  can  be  either  resolved  or  parameterized  via  detailed 
measurement  and  modeling.  These  new  initiatives  also  have  a  similar  goal  to  establish  improved 
understanding  on  coastal  hydrodynamics  and  sediment  transport  processes  in  heterogeneous 
environment  with  more  emphasis  on  cohesive  sediment  transport. 

The  dynamics  of  fluid-mud  transport  involve  a  variety  of  physical  mechanisms,  including  for  example, 
the  boundary  layer  and  gravity-driven  transport,  turbulence  modulation,  flocculation,  non-Newtonian 
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rheological  behavior  and  consolidation  (Dyer  1989;  Mehta  1989).  Hence,  a  general  modeling 
framework  for  fluid  mud  needs  to  be  based  on  multiphase  flow  theory.  In  this  study,  a  fine  sediment 
transprt  modeling  framework  based  on  Fast  Equilibrium  Eulerian  Approximation  (Ferry  and 
Balachandar  2001)  to  the  multiphase  equations  has  been  developed  and  extended  to  model  various 
cohesive  sediment  transport  processes.  The  cohesive  sediment  modeling  framework  is  the  basis  of  the 
three  numerical  models  developed  for  IDV,  2DV  and  full  3D  for  different  applications.  In  the  IDV 
modeling,  the  dynamics  of  wave-supported  gravity-driven  mudflows  is  studied  (Hsu  et  al.  2007;  2008). 
In  the  2DV  formulation,  wave-mud  interaction  and  sedimentation  at  salinity  stratified  river  mouth  are 
simulated.  To  resolve  detailed  turbulence-sediment  interactions  and  to  improve  the  closures  in  the 
Reynolds-averaged  approach.  Direct  Numerical  Simulation  (DNS)  of  fluid  mud  transport  under 
oscillatory  flow  is  carried  out  in  a  full  3D  simulatoin  based  on  pseudo-spectrum  scheme. 

WORK  COMPLETED 

Wave-supported  gravity-driven  mudflows:  Typical  coastal  models,  such  as  the  ongoing  NOPP- 
Community  Sediment  Transport  Model  (NOPP-CSTM,  Warner  et  al.  2008)  are  not  designed  to  resolve 
the  thin  wave  boundary  layer  near  the  seabed.  Hence,  processes  that  occur  within  the  wave  boundary 
layer  need  to  be  parameterized  as  sub-module  to  provide  quantities  such  as  bottom  drag  coefficient  and 
sediment  transport  rate.  The  IDV  fluid  mud  model  is  utilized  to  study  fluid  mud  transport  in  the  wave- 
current  boundary  layer.  The  characteristics  of  wave-supported  gravity-driven  mudflows  are  diagnosed 
by  varying  the  bed  erodibility,  floe  properties  (fractal  dimension),  rheological  stresses,  and  floe 
dynamics  in  the  numerical  simulations.  Model  results  for  moderate  concentration  suggest  that  using  an 
appropriately  specified  fractal  dimension,  the  dynamics  of  wave-supported  gravity-driven  mudflow 
can  be  predicted  without  explicitly  incorporating  rheological  stress.  However,  incorporating 
rheological  stress  makes  the  results  less  sensitive  to  prescribed  fractal  dimension.  For  high 
concentration  condition,  it  is  necessary  to  incorporate  rheological  stress  in  order  to  match  observed 
intensity  of  downslope  gravity-driven  current.  Model  results  are  further  analyzed  to  evaluate  and 
calibrate  simple  parameterizations.  Analysis  suggests  that  when  neglecting  rheological  stress,  the  drag 
coefficient  decreases  with  increasing  wave  intensity,  and  seems  to  follow  a  power  law.  However,  when 
rheological  stress  is  incorporated,  the  resulting  drag  coefficient  is  more  or  less  constant  (around 
0.0013)  for  different  wave  intensities.  Model  results  further  suggest  the  bulk  Richardson  numbers  has  a 
magnitude  smaller  than  0.25  and  is  essentially  determined  by  the  amount  of  available  soft  mud  (i.e., 
the  erodibility),  suggesting  a  supply-limited  condition  for  unconsolidated  mud.  In  order  to  predict  more 
quantitatively  the  dynamics  of  wave-supported  gravity-driven  mudflows,  more  detailed  information  on 
floe  properties,  rheological  stress  and  erodibility  are  necessary.  Findings  of  this  study  are  summarized 
in  a  manuscript  submitted  for  publication. 

Flocculation;  The  present  cohesive  sediment  transport  framework  requires  realistic  constitutive 
equation  for  floe  breakup  and  aggregation  processes  (e.g.,  Winterwerp  1998).  Based  on 
field/laboratory  observations,  the  fractal  dimension  is  not  a  constant  and  shall  change  dynamically  with 
the  carrier  turbulent  flow  (e.g.,  Dyer  and  Manning  1999).  Recently,  Khelifa  and  Hill  (2006)  proposed  a 
model  for  settling  velocity  that  utilizes  a  variable  fractal  dimension  depending  on  the  floe  size  itself 
(e.g.,  fractal  dimension  becomes  3.0  when  floe  size  reduces  to  that  of  primary  particle).  Here,  we 
further  derive  the  balance  equations  for  floe  dynamics  based  on  the  variable  fractal  dimension  and 
come  up  a  new  model  for  floe  dynamics  that  can  be  coupled  with  numerical  model  for  fluid  mud 
transport.  The  new  formulation  is  shown  to  predict  the  equilibrium  floe  size  under  different  turbulent 
shear  rates  measured  in  the  laboratory  (Son  and  Hsu  2008).  Currently,  we  are  testing  and  accuracy  and 
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robustness  of  such  flocculation  formulations  in  a  numerical  model  for  cohesive  sediment  transport  in 
tidal-dominated  condition. 

Wave-mud  interaction:  As  a  part  of  this  framework,  a  well- validated  2DV  depth/phase-resolving 
wave  model  (COBRAS,  Lin  &  Liu,  1998)  based  on  the  Reynolds-averaged  Navier-Stokes  (RANS) 
equations  has  been  modified  for  the  study  of  wave-mud  interaction.  The  resulting  governing  equations 
and  closures  reduce  to  the  RANS  equations  when  the  sediment  concentration  approaches  zero.  Hence, 
the  numerical  model  is  able  to  simulate  continuously  and  consistently  the  nonlinear  water  wave 
propagation,  the  fluid-mud  generation  and  transport,  the  wave-boundary  layer  processes,  turbulence 
modulation  owing  to  the  presence  of  the  fluid-mud,  and  the  rheological  effects  on  attenuating  the 
waves  with  a  single  set  of  balance  equations  and  closures.  Since  FY08,  postdoc  researcher  Alec 
Torres-Freyermuth  has  worked  closely  with  Hsu  and  Traykovski  as  well  as  other  teams  participated  in 
the  Achafalaya  field  exepriment  to  utilize  this  2DV  model  to  study  wave-mud  interaction.  The  code 
has  been  further  improved  upon  the  implementation  of  a  new  wave  forcing  boundary  condition 
(Torres-Freyermuth  et  al.  2007),  the  optimization  of  the  mesh  and  obstacle  generation,  and  the  CPU 
time  reduction  that  allow  computing  more  realistic  and  longer  simulations. 

RESULTS 

With  the  aim  of  gaining  insight  on  the  processes  related  to  wave-mud  interaction,  the  2DV  model  is 
employed  using  three  different  wave  types.  Cohesive  sediment  of  floe  size  22  pm  and  specific  gravity 
1.34  is  specified  in  all  tests  shown  here.  High  resolution  in  the  vertical  water  column  is  adopted  and 
hence  wave  boundary  layer  structure  and  luctocline  are  resolved  by  the  numerical  model. 
Monochromatic  progressive  waves,  standing  waves,  and  wave  groups  propagating  over  a  planar 
bottom  are  considered.  The  monochromatic  progressive  wave  cases  are  studied  as  a  part  of  the  pre¬ 
validation  stage  of  the  model.  On  the  other  hand,  the  standing-wave  and  wave  group  cases  are  further 
studied  in  order  to  test  different  hypothesis.  The  former  concerns  the  scaling  between  the  wave  bottom 
boundary  layer  and  the  fluid  mud  layer  thickness  reported  in  the  field  (e.g.  Traykovski  et  al.  2007), 
while  the  later  is  used  to  obtain  more  insight  on  the  role  of  nonlinear  wave  interaction  in  the  damping 
of  the  higher  (e.g.  Sheremet  &  Stone,  2003)  and  the  lower  frequencies  waves  (e.g.  Elgar  & 
Raubenheimer,  2008)  in  the  incident  wave  spectrum. 

For  the  monochromatic  wave  tests,  a  cnoidal  wave  train  of  0.72  m  wave  height  and  6  s  period  is  sent 
into  a  2.5  m  water  depth  numerical  flume.  The  amount  of  suspended  sediment  is  controlled  by  varying 
the  sediment  availability  (erodibility)  for  each  of  the  two  cases,  resulting  in  different  fluid-mud 
characteristics  (e.g.  thickness,  viscosity,  etc.)  and  damping  rates  under  the  same  forcing.  Subsequently, 
the  numerically  obtained  time-  and  depth-dependent  fluid-mud  layer  characteristics  are  integrated  to 
provide  the  input  bulk  parameters  (mud  layer  thickness  5  and  bulk  viscosity  Um)  required  in  simplified 
two-layer  wave  propagation  model  (Kranenburg,  2008).  Good  agreement  between  the  2DV  model 
results  and  the  two-layer  model  results  is  observed  in  terms  of  the  wave  height  at  four  different  cross¬ 
shore  positions  (Figure  lb).  These  results  provide  confidence  in  terms  of  the  model  consistency  and 
capability  for  providing  reliable  information  on  the  magnitude  of  wave  dissipation  rates  and  the  model 
application  on  the  study  of  more  complicated  wave  conditions  (e.g.  wave  groups,  irregular  waves). 

A  standing  wave  test  is  selected  to  study  bottom  wave  boundary  layer  dynamics.  A  small  amplitude 
monochromatic  long-wave  is  propagated  into  the  numerical  wave  flume  and  is  completely  reflected  by 
a  vertical  wall  located  at  the  downstream  end  of  the  flume  generating  a  standing  wave  pattern  (Figure 
2a).  Fluid-mud  is  preferentially  suspended  near  the  nodes  (horizontal  velocity  amplitude  maximum) 
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and  in  fact  the  fluid-mud  thickness  modulates  in  the  cross-shore  direction  depending  on  the  horizontal 
velocity  amplitude  (Figure  2b).  The  wave  boundary  layer  thickness  (blue  curve  in  Figure  2b)  coincides 
with  the  luctocline  location  defined  as  the  largest  vertical  gradient  in  the  sediment  concentration  profile 
(red  curves  in  Figure  2b).  Moreover,  the  velocity  profiles  with  fluid-mud  (blue  curves)  when  compared 
with  the  clear  fluid  condition  (black-dashed  curve)  shows  a  significant  wave  boundary  layer 
enhancement.  In  terms  of  the  wave  amplitude,  comparison  between  the  clear  fluid  and  the  fluid-mud 
results  shows  evidence  of  wave  energy  transfer  from  principle  harmonic  to  higher  harmonics  due  to 
wave-wave  interaction  triggered  by  the  presence  of  fluid-mud  (Figure  2c).  The  reader  is  referred  to 
unexpected  local  amplitude  maximum  near  the  node  when  the  fluid-mud  is  present. 

Finally,  a  wave  group  with  its  corresponding  sub-  and  super-harmonics  is  sent  to  propagate  into  the 
numerical  wave  flume.  A  sponge  layer  is  placed  at  the  downstream  end  of  the  flume  in  order  to  absorb 
the  incoming  wave  energy.  Free-surface  elevation  time  series  measured  downstream  of  the  flume 
shows  a  wave  decay  as  compared  with  another  simulation  without  the  fluid  mud  (Figure  3a).  The 
energy  attenuation  occurs  almost  among  all  frequencies  (except  around  f=0.1  &  0.35  Hz)  in  the  energy 
spectrum  (Figure  3b).  The  frequency  dependence  of  the  dissipation  rate  is  obtained  by  comparing  the 
numerical  simulation  from  the  clear  fluid  with  the  one  with  fluid-mud.  Differences  between 
nondissipative  (blue)  and  dissipative  (red)  model  predictions  are  attributed  to  mud-induced  dissipation 
and  mud-induced  wave  energy  transfer.  The  higher  dissipation  rates  for  this  case  correspond  to  the 
sub-  and  super-harmonic  frequencies  (Figure  3c)  consistent  with  previous  observations  pointing  out  the 
importance  of  nonlinear  energy  transfer  as  a  mechanism  of  wave  attenuation  (e.g.  Elgar  & 
Raubenheimer  2008).  Moreover,  the  energy  increase  at  discrete  frequencies  (e.g.  f=0. 1  Hz  in  Figure 
3b)  for  the  fluid-mud  simulation  is  another  indication  of  the  excitation  of  wave-wave  interaction  by  the 
presence  of  mud. 

IMPACT/APPLICATIONS 

The  present  research  efforts  focus  on  developing  a  numerical  modeling  framework  for  cohesive 
sediment  transport.  Currently,  the  Pis  are  also  actively  participating  in  other  related  initiatives  and  are 
partners  of  the  ongoing  NOPP-CSTM  and  CSDMS.  Hsu  has  actively  participated  in  the  ongoing 
NICOP  -  Preliminary  Investigations  on  the  Fate  of  Terrestrial  Sediments  in  the  Coastal  Ocean 
Discharged  from  Taiwanese  Small  Mountainous  Rivers,  supported  by  ONRG  and  ONR  (Coastal 
Geosciences).  The  2DV  model  discussed  here  was  recently  extended  to  calculate  salinity  transport  to 
simulate  sedimentation  processes  at  river  mouth. 


4 


20  40  60  80  100  120  140  160  180 


x[m] 


Figure  1:  (a)  Snapshot  of  the  model  results  for  free-surface  elevation  and  mud  concentration  in  the 
numerical  wave  flume.  The  vertical  red-dashed  lines  correspond  to  the  discrete  cross-shore  positions 
where  wave  gauge  data  was  collected,  (b)  Comparison  of  wave  gauge  data  (symbols)  and  the  two- 
layer  model  predictions  (lines)  for  wave  height  as  a  function  of  cross-shore  position.  The  triangles 
(circles)  and  red-dashed  (blue-solid)  line  represent  2DV model  results  and  two-layer  model  results 
respectively  for  the  low  (high)  erodibility  condition  with  bulk  mud  layer  characteristic  required  by 
the  two-layer  model  specified  as  d=6cm,  Vm=7.5xl(T‘^m^/s,  and  pm=1300Kg/m^  (d=12cm,  Vm=3.8xlff' 
m  /s,  and  pm=1300Kg/m  ).  These  bulk  mud  layer  parameters  are  obtained  directly  by  integrating 

2DV model  results. 
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Figure  2:  (a)  Snapshot  of  the  model  results  for  free-surface  elevation  and  mud  concentration  in  the 
numerical  wave  flume  for  the  standing  wave  case.  The  dashed  square  box  represents  the  area 
discussed  in  (b)  and  (c).  (b)  Snapshot  of  fluid-mud  concentration  near  the  bed  with  the 
corresponding  velocity  (blue  curve)  and  concentration  (red  curve)  profiles  under  the  antinode  (x=21 
m)  and  the  node  (x=34  m).  The  dotted  curves  represent  zero  mud  concentration  which  is  a  reference 
case  of  clear  fluid  simulated  for  the  sake  of  comparing  the  effect  of  mud  with  black-dashed  curves 
represent  the  corresponding  velocity  profiles.  The  wave  boundary  layer  thickness  scales  with  the 
fluid-mud  layer,  and  is  enhanced  when  comparing  to  the  clear  fluid  condition,  (c)  Comparison  of 
the  cross-shore  variation  of  wave  amplitude  for  the  clear  fluid  (blue)  and  the  fluid-mud  (red)  tests. 

The  first  harmonic  (solid  curves)  is  attenuated  due  to  the  presence  of  the  mud  at  all  cross-shore 
positions.  On  the  other  hand,  the  second-  (dashed  curves)  and  third-  harmonics  (dotted  curves)  are 
enhanced  neat  the  node  (x=34  m)  when  fluid  mud  presents.  Hence  the  total  wave  amplitude  (sum  of 
all  curves)  near  the  node  (x=34m)  is  unexpectedly  enhanced. 
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Figure  3:  (a)  Free-surface  elevation  time  series  measured  at  the  downstream  of  the  flume  (x=185m) 
for  the  clear  fluid  (blue  curve)  and  the  fluid-mud  (red  curve)  tests,  (b)  Free-surface  energy  density 
spectra  for  these  two  simulations.  Differences  between  the  non-dissipative  (blue  curve)  and 
dissipative  (red  curve)  cases  are  attributed  to  fluid-mud  induced  dissipation  and  mud  induced  wave 
energy  transfer,  (c)  Dissipation  rate  at  discrete  frequencies  obtained  as  the  difference  in  the  energy 
flux  from  the  two  simulations.  The  dissipation  rates  are  higher  for  the  sub-harmonic  (f2-fl)  and 

super-harmonic  (f2+fl)  frequencies. 
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